您好,欢迎来到达州工业云! 平台首页 企业驾驶舱 帮助中心 企业登录 企业注册

HI,欢迎使用达州工业云平台!

账号必须大于2位

创新资源平台
服务平台首页>专利库>专利详情

ISAR图像距离单元选择横向定标方法

  • 申请号:CN201210476554.X 申请公布号: CN102944875B
  • 申请日: 2012-11-21 申请公布日: 2014-07-23
  • 申请(专利权)人: 专利代理机构: 陕西电子工业专利中心
  • 分类号:G01S13/89

专利介绍

本发明公开一种ISAR图像距离单元选择横向定标方法,主要解决现有方法没有选择到合适的距离单元以至于调频率估计误差较大的问题。其实现过程是:读取雷达脉冲回波;对雷达回波数据进行运动补偿;根据幅度选取目标所在的距离单元;从目标所在的距离单元中根据多种准则选择所需距离单元;利用局部多项式傅里叶变换对选择的距离单元进行调频率估计;根据归一化面积准则二次选择距离单元;多次利用最小二乘法求平均,得到角速度估计值,由角速度完成对ISAR图像的横向定标。本发明提高了调频率估计的准确性,可用于对ISAR图像的横向定标。
1.一种ISAR图像距离单元选择横向定标方法,包括如下步骤:(1)获取雷达脉冲回波信号为SR(m,n),并对该雷达脉冲回波信号进行运动补偿,得到运动补偿后的雷达脉冲回波信号,表示为ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="200"/>其中,m为快时间离散采样,m=1,2,...,M,M为距离单元数,n为慢时间离散采样,n=1,2,...,N,N为脉冲数;(2)根据雷达脉冲回波信号ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="175"/>的幅度信息选取ISAR图像中目标所在的距离单元,得到目标回波信号ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="186"/>其中,l=i,i+1,...,j,l表示目标所在的距离单元,i表示目标所在距离单元的起始位置,j表示目标所在距离单元的末尾位置;(3)利用多种准则从雷达脉冲回波信号ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="172"/>中选择所需的距离单元:(3a)利用距离-多普勒算法对目标回波信号ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="151"/>成像,得到ISAR图像I(l,n"),其中n"表示频率,n"=1,2,...,N,分别计算ISAR图像I(l,n")每个距离单元的对比度c(l): c ( l ) = mean { [ | I ( l , n " ) | - mean [ | I ( l , n " ) | ] ] 2 } mean [ | I ( l , n " ) | ] , ]]>其中,|·|表示求模运算;(3b)将对比度c(l)进行降序排序,根据设定的对比度阈值μ1选择出c(l)的前L个值,从目标回波信号ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="155"/>中提取这L个距离单元,作为初选距离单元,记为S"(l",n),其中,l"是初选出来的L个距离单元,l"=1,2,...,L;(3c)计算初选距离单元S"(l",n)的每个距离单元的对比度c(l");(3d)分别计算初选距离单元S"(l",n)每个距离单元的包络归一化方差v(l"),相位匹配二项式误差p(l")以及每个距离单元的平均幅度a(l"): v ( l " ) = cov ( S " ( l " , n ) ) cov ( S " ( l " , n ) ) + [ mean ( S " ( l " , n ) ) ] 2 ]]> p ( l " ) = Σ n | angle ( S " ( l " , n ) ) - LSE ( angle ( S " ( l " , n ) ) ) | ]]> a ( l " ) = Σ n ( S " ( l " , n ) ) N , ]]>其中,cov(·)表示求方差运算,mean(·)表示求平均值运算,angle(·)表示求相位运算,LSE(·)表示最小二乘拟合直线,∑(·)表示求和运算;(3e)分别对包络归一化方差v(l")进行升序排序,对相位匹配二项式误差p(l")进行升序排序,对每个距离单元的平均幅度a(l")进行降序排序,对每个距离单元的对比度c(l")进行降序排序,依据排序后的v(l"),p(l"),a(l")和c(l")分别对初选距离单元S"(l",n)中每个距离单元进行“投票”:(3e1)找到初选距离单元S"(l",n)中每个距离单元在v(l"),p(l"),a(l")和c(l")中的排名顺序,分别对每个距离单元在v(l"),p(l"),a(l")和c(l")中排名顺序的名次进行相加,计算得到每个距离单元的名次总和;(3e2)对每个距离单元的名次总和重新进行升序排序,选择出前L"个值,从初选距离单元S"(l",n)中选择出这L"个值对应的L"个距离单元,作为“投票”选择出的距离单元,记为S"(l"",n),其中,l""是新选择出来的L"个距离单元,l""=1,2,...,L",L" IC ( l " " , α ) = mean { [ lpft ( l " " , f , α ) - mean [ lpft ( l " " , f , α ) ] ] 2 } mean [ lpft ( l " " , f , α ) ] ]]>其中,lpft(l"",f,α)为局部多项式傅里叶变换,计算公式为 lpft ( l " " , f , α ) = | Σ n = 0 fn S " ( l " " , n ) e - j 2 π ( f ( n T r ) + α 2 ( n T r ) 2 ) | , ]]>该式可采用快速傅里叶变换实现,f表示频率,Tr表示脉冲重复周期,mean[lpft(l"",f,α)]表示在频率维f上的求平均值运算;(4c)根据图像对比度IC(l"",α)计算调频率chirp(l"")和每个调频率对应的归一化面积area(l""), chirp ( l " " ) = arg max α [ IC ( l " " , α ) ] area ( l " " ) = Σ α ( IC ( l " " , α ) ) , ]]>其中,ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="215"/>表示对表达式求最大时对应的α;(4d)更新l""=l""+1,判断l""=N是否成立;若成立则继续进行,否则转到步骤(4b);(5)对归一化面积area(l"")进行升序排序,根据设定的归一化面积阈值μ2选择出area(l"")的前K个值,从调频率chirp(l"")中得到这K个值所对应的调频率chirp(k),从“投票”选择出的距离单元ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="159"/>中得到其对应的K个距离单元,从雷达回波信号ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="174"/>中得到这K个距离单元对应的距离单元位置rangebin(k),其中,k表示最后选择出来的距离单元,k=1,2,...,K;(6)根据最小二乘方法,由前K个调频率chirp(k)和其对应的距离单元rangebin(k)线性拟合斜率ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="88"/>由前K-1个调频率chirp(k-1)和其对应的距离单元rangebin(k-1)线性拟合斜率ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="119"/>依次类推,直到由前2个调频率chirp(2)和其对应的距离单元rangebin(2)线性拟合,得到斜率ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="81"/>最后计算所求斜率的平均值ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="64"/>即 a ^ = 1 K - 1 Σ k = 2 K a ^ k ; ]]>(7)根据斜率平均值ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="66"/>得到垂直于成像平面的角速度ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="78"/> Ω ^ = c 2 f 0 | a ^ | ]]>其中,c为光速,f0为雷达射频;(8)依据角速度ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="66"/>对步骤(3a)中得到的ISAR图像I(l,n")的方位维坐标n"乘上方位分辨率ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="159"/>即可完成ISAR图像的横向定标;再对ISAR图像I(l,n")的距离维坐标l乘上距离分辨率ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="102"/>即可得到定标后的ISAR图像,表示为:ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="429"/>其中,B为雷达信号带宽,Fr为脉冲重复频率,λ为波长。
2.根据权利要求1所述的一种ISAR图像距离单元选择横向定标方法,其中步骤(6)所述的根据最小二乘方法,由前K个调频率chirp(k)和其对应的距离单元rangebin(k)线性拟合,得到斜率ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="86"/>通过如下计算公式进行: a ^ k = Σ k = 1 K chirp ( k ) rangebin ( k ) - 1 K Σ k = 1 K chirp ( k ) Σ k = 1 K rangebin ( k ) Σ k = 1 K ( rangebin ( k ) ) 2 - 1 K ( Σ k = 1 K rangebin ( k ) ) 2 . ]]>技术领域
本发明属于雷达技术领域,涉及ISAR横向定标方法,特别是利用多种准则选择距离单元,可用在基于ISAR图像的目标识别,对于其中的特征提取具有重要的支撑作用。
背景技术
ISAR横向定标是基于ISAR成像的目标识别分类中涉及的关键技术。横向定标算法的核心就在于获得在成像时间内目标相对雷达视线转过的总转角,所有的横向定标的算法也是围绕这一点展开。目前横向定标算法大致可分为以下两类:一是基于图像域的定标算法,包括图像熵法,基于对称目标的横向定标算法,基于图像旋转自相关的横向定标算法以及基于2Dxa0FT与极化映射的横向定标算法等。其中图像熵法采用评价函数来优化未知的转角参数,当图像评价函数获得极值,即图像最聚焦时,得到相应的转角,但这种方法需要在解空间内反复搜索,计算量非常大;基于对称目标的横向定标算法是基于对称目标的对称轴垂直于对称点直线,利用主轴直线与对称点直线斜率之间的关系得到目标的转角估计,但是该方法只能针对对称目标,而且主轴直线与对称点直线斜率的准确获取存在一定困难;基于图像旋转的横向定标算法利用图像旋转技术将一幅图像经过旋转后与另一幅图像进行匹配,两幅图像的匹配程度可以通过两幅图像的相关系数来衡量,相关系数越大,图像的匹配程度就越高。在匹配程度最高时得到对应的转角,即可实现对目标ISAR图像的横向定标,但该算法需要较高质量的ISAR图像,而且需要确定图像的旋转中心;基于2Dxa0FT与极化映射的横向定标算法则对图像旋转中心问题进行了优化,即利用二维傅里叶变换将图像变换到频域,使该频域图像的旋转中心处于图像的几何中心,然后利用频域图像的旋转相关性得到目标转角,从而实现横向定标,但该算法同样对ISAR图像质量要求较高,而且计算量也较大。二是基于雷达回波数据域的定标算法,包括轨迹拟合法,特显点跟踪法,最大似然参数估计法,以及基于LPFT的调频率估计。其中轨迹拟合法主要是利用雷达对目标跟踪的窄带数据,得到目标的运动轨迹,然后得到运动目标相对雷达的转角变化,但由于窄带数据的误差较大,以及目标的翻滚,倾斜和偏航等复杂运动,导致这种方法的估计精度非常低;特显点跟踪法是分别跟踪目标上3个不同的孤立散射点,但实际情况下很难找到高质量的孤立点;最大似然参数估计法是通过最大似然函数来估计回波信号的三次相位系数和一次相位系数,然后利用两者之间的关系得到目标的转角估计,但由于三次相位系数一般非常小,且由于同一距离单元中散射点之间的相互影响,导致实际情况下很难准确估计回波信号中的三次相位系数,而且该方法涉及多维搜索,计算量也较大;基于LPFT的调频率估计则是利用LPFT变换搜索回波信号的二次调频率,根据二次调频率与目标转速之间的关系得到目标在成像时间内的总转角,但该算法对距离单元的选择较为敏感。基于LPFT的调频率估计是目前较为良好的横向距离定标算法。但是它的不足之处也较为明显,如果调频率估计不准,则通过最小二乘拟合难以得到正确的斜率,角速度的估计也随之错误,从而造成横向距离定标的结果误差较大。在实测数据中,由于噪声、目标运动状态、相邻特显点的相互作用等各种因素的存在,对调频率估计的结果影响很大,事实上通过很多距离单元估计的调频率是错误的。如何选择合适的距离单元进行调频率估计是此算法在实际应用中最为关键的问题。
发明内容
本发明的目的在于上述已有技术的不足,提出一种ISAR图像距离单元选择横向定标方法,以提高调频率估计的准确性。本发明的技术方案是:以基于调频率估计的ISAR横向距离定标方法为基础,设计多种准则选择距离单元和对估计的调频率进行二次筛选,其实现步骤包括如下:(1)获取雷达脉冲回波信号为SR(m,n),并对该雷达脉冲回波信号进行运动补偿,得到运动补偿后的雷达脉冲回波信号,表示为S′R(m,n),其中,m为快时间离散采样,m=1,2,...,M,M为距离单元数,n为慢时间离散采样,n=1,2,...,N,N为脉冲数;(2)根据雷达脉冲回波信号S′R(m,n)的幅度信息选取ISAR图像中目标所在的距离单元,得到目标回波信号S′R(l,n),其中,l=i,i+1,...,j,l表示目标所在的距离单元,i表示目标所在距离单元的起始位置,j表示目标所在距离单元的末尾位置;(3)利用多种准则从雷达脉冲回波信号S′R(m,n)中选择所需的距离单元。(3a)利用距离-多普勒算法对目标回波信号S′R(l,n)成像,得到ISAR图像I(l,n′),其中n′表示频率,n′=1,2,...,N,分别计算ISAR图像I(l,n′)每个距离单元的对比度c(l): c ( l ) = mean { [ | I ( l , n ′ ) | - mean [ | I ( l , n ′ ) | ] ] 2 } mean [ I ( l , n ′ ) ] , ]]>其中,|·|表示求模运算;(3b)将对比度c(l)进行降序排序,根据设定的对比度阈值μ1选择出c(l)的前L个值,从目标回波信号S′R(l,n)中提取这L个距离单元,作为初选距离单元,记为S′(l′,n),其中,l′是初选出来的L个距离单元,l′=1,2,...,L;(3c)计算初选距离单元S′(l′,n)的每个距离单元的对比度c(l′);(3d)分别计算初选距离单元S′(l′,n)每个距离单元的包络归一化方差v(l′),相位匹配二项式误差p(l′)以及每个距离单元的平均幅度a(l′): v ( l ′ ) = cov ( S ′ ( l ′ , n ) ) cov ( S ′ ( l ′ , n ) ) + [ mean ( S ′ ( l ′ , n ) ) ] 2 ]]> p ( l ′ ) = Σ n | angle ( S ′ ( l ′ , n ) ) - LSE ( angle ( S ′ ( l ′ , n ) ) ) | ]]> a ( l ′ ) = Σ n ( S ′ ( l ′ , n ) ) N , ]]>其中,cov(·)表示求方差运算,mean(·)表示求平均值运算,angle(·)表示求相位运算,LSE(·)表示最小二乘拟合直线,∑(·)表示求和运算;(3e)分别对包络归一化方差v(l′)进行升序排序,对相位匹配二项式误差p(l′)进行升序排序,对每个距离单元的平均幅度a(l′)进行降序排序,对每个距离单元的对比度c(l′)进行降序排序,依据排序后的v(l′),p(l′),a(l′)和c(l′)分别对S′(l′,n)中每个距离单元进行“投票”,由此提取出S′(l′,n)中前L′个距离单元,作为“投票”选择出的距离单元,记为S′(l′′,n),其中,l′′是新选择出来的L′个距离单元,l′′=1,2,...,L′,L′<L。(4)对“投票”选择出的距离单元S′(l′′,n)的每个距离单元进行局部多项式傅里叶变换lpft,得到调频率chirp(l′′)和每个调频率对应的归一化面积area(l′′):(4a)设定参数:α=-1,-0.999,-0.998,...,1,fn=16*N;其中,α表示调频率的搜索范围,fn表示傅里叶变换的点数;(4b)对“投票”选择出的距离单元S′(l′′,n)中第l′′个距离单元进行局部多项式傅里叶变换lpft(l′′,f,α),并计算第l′′个距离单元的图像对比度IC(l′′,α),并令l′′=1, IC ( l ′ ′ , α ) = mean { [ lpft ( l ′ ′ , f , α ) - mean [ lpft ( l ′ ′ , f , α ) ] ] 2 } mean [ lpft ( l ′ ′ , f , α ) ] , ]]>其中,lpft(l′′,f,α)为局部多项式傅里叶变换,计算公式为ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="700"/>该式可采用快速傅里叶变换实现,f表示频率,Tr表示脉冲重复周期,mean[lpft(l′′,f,α)]表示在频率维f上的求平均值运算;(4c)根据图像对比度IC(l′′,α)计算调频率chirp(l′′)和每个调频率对应的归一化面积area(l′′), chirp ( l ′ ′ ) = arg max α [ IC ( l ′ ′ , α ) ] ]]> area ( l ′ ′ ) = Σ α ( IC ( l ′ ′ , α ) ) , ]]>其中,ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="195"/>表示对表达式求最大时对应的α;(4d)更新l′′=l′′+1,转到步骤(4b),直至l′′=N。(5)对归一化面积area(l′′)进行升序排序,根据设定的归一化面积阈值μ2选择出area(l′′)的前K个值,从调频率chirp(l′′)中得到这K个值所对应的调频率chirp(k),从“投票”选择出的距离单元S′(l′′,n)中得到其对应的K个距离单元,从雷达回波信号S′R(m,n)中得到这K个距离单元对应的距离单元位置rangebin(k),其中,k表示最后选择出来的距离单元,k=1,2,...,K;(6)根据最小二乘方法,对调频率chirp(k)和其对应的距离单元rangebin(k)线性拟合,得到K个调频率估计的斜率ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="63"/>对调频率chirp(k-1)和其对应的距离单元rangebin(k-1)线性拟合,得到K-1个调频率估计的斜率ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="89"/>依次类推,直到对调频率chirp(2)和其对应的距离单元rangebin(2)线性拟合,得到2个调频率估计的斜率ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="59"/>最后计算所求斜率的平均值ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="46"/>即ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="314"/>(7)根据斜率平均值ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="44"/>得到垂直于成像平面的角速度ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="57"/> Ω ^ = c 2 f 0 | a ^ | ]]>其中,c为光速,f0为雷达射频;(8)依据角速度ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="57"/>对步骤(3a)中得到的ISAR图像I(l,n′)的方位维坐标n′乘上方位分辨率ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="136"/>即可完成ISAR图像的横向定标;再对ISAR图像I(l,n′)的距离维坐标l乘上距离分辨率ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="91"/>即可得到定标后的ISAR图像,表示为:ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="410"/>其中,B为雷达信号带宽,Fr为脉冲重复频率,λ为波长。本发明具有以下优点:1)本发明基于雷达回波数据域进行ISAR图像横向距离定标,相比传统基于图像域的定标方法更加稳健;2)本发明通过多种准则选择距离单元,最大程度的去掉不适合做调频率估计的距离单元,使得调频率估计结果更加精确;3)本发明通过对不同数目的距离单元多次估计的调频率进行平均,提高了调频率估计的稳健性。
附图说明
图1是本发明的流程图;图2是本发明对实测数据的成像结果图;图3是本发明对目标所在距离单元进行初选的结果图;图4是用本发明方法选择的距离单元以及调频率的最小二乘拟合图;图5是用常规方法选择的距离单元以及调频率的最小二乘拟合图。
具体实施方式
参照图1,本发明的具体实现步骤如下:步骤1,将雷达脉冲回波表示为SR(m,n),其中蒯为快时间离散采样,m=1,2,...,M,M为距离单元数;n为慢时间离散采样,n=1,2,...,N,N为脉冲数。步骤2,对雷达脉冲回波进行运动补偿,得到运动补偿后的信号,表示为S′R(m,n)。(2a)采用相邻相关法进行包络对齐。(2a1)将第n次和第n+1次相邻两次回波分别表示为SR(m,n)和SR(m,n+1),它们的互相关函数为:Rn,n+1(k)=∑SR(m,n)SR(m-k,n+1);利用该互相关函数Rn,n+1(k),计算Rn,n+1(k)峰值所对应的时延值k;在实际计算互相关函数时,通常将时间离散值m做8倍的插值处理;(2a2)利用时延值k对SR(m,n+1)进行包络平移;(2b)对包络对齐后的数据进行相位校正。(2b1)使用多普勒中心法将第n和第n+1次相邻回波信号共轭相乘,按幅度乘积进行加权平均,得到其相位差指数函数,即按下式求出多普勒中心相位差的复指数函数ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="96"/>ntent="drawing" img-format="tif" inline="no" orientation="portrait" wi="417"/>其中,Sm,n,Smn+1表示第n和第n+1次相邻回波,ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="70"/>表示Sm,n的共轭;(2b2)把得到的复指数函数ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="74"/>补偿到每一距离单元上,即可实现初步的相位校正,初步校正后的数据记为ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="96"/>(2b3)对初步校正的数据,继续使用最小熵单调迭代算法进行相位校正,通过迭代计算出相邻回波的相位误差φ,表示为φ={φ1,φ2,...,φn,...,φN},φn表示补偿到距离单元第n个位置的相位因子,n=1,2,...,N,并将第a次迭代估计的相位误差记为φ(a),其中ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="651"/>a为迭代次数,ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="64"/>表示第a次补偿到距离单元第n个位置的相位因子,得到复指数函数ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="109"/> e j φ ( a ) = { e j φ 1 ( a ) , e j φ 2 ( a ) , . . . , e j φ n ( a ) , . . . , e j φ N ( a ) } ; ]]>(2b4)把复指数函数ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="81"/>补偿到每一距离单元上,即可实现数据的相位校正,具体步骤如下:i)初始化a=0,φ(0)={0,0,...,0};i)利用相位误差估计值φ(a)对初步相位校正的数据ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="73"/>进行相位补偿,得到初步聚焦处理的ISAR图像ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="92"/> a m , k ( a ) = 1 / N Σ n = 0 N - 1 y ~ m , n e j φ n ( a ) e - j 2 πkn / N , ]]>其中,k表示频率维的离散采样,k=1,2,...,N;iii)对第a次补偿到距离单元第n个位置的相位因子ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="64"/>进行迭代求解,得到第a+1次补偿到距离单元第n个位置的相位因子ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="118"/>更新第a次迭代估计的相位误差φ(a),即: φ n ( a + 1 ) = φ n ( a ) + tan - 1 ( B n / A n ) , ]]>其中,An,Bn是中间变量, A n = 2 NE Z Im { Σ m y ~ m , n Σ k e j 2 πkn / N z m , k * ln | z m , k ( a ) | 2 } , ]]> B n = 2 NE Z Re { Σ m y ~ m , n Σ k e j 2 πkn / N z m , k * ln | z m , k ( a ) | 2 } - 2 N 2 E Z Σ m { | y ~ m , n | 2 Σ k ln | z m , k ( a ) | 2 } , ]]>式中,Im(·)是求复数的虚部运算,Re(·)是求复数的实部运算,ln(·)表示求自然对数运算,∑(·)表示求和运算,Ez为图像能量,ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="684"/>ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="65"/>为zm,k的共轭,zm,k是重新聚焦处理后的ISAR图像,ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="602"/> φ n = φ n ( a ) ; ]]>iv)更新a=a+1,转到步骤ii),直到图像能量Ez满足要求。本步骤使用的最小熵单调迭代算法基于ISAR图像熵函数,寻找ISAR图像熵函数的替代函数,优化求解过程,其详细介绍见文献“Thomasxa0J.Kragh,′′MonotonicIterativexa0Algorithmxa0Forxa0Minimum-Entropyxa0Autofocus",Lincolnxa0Laboratory,MIT”。步骤3,根据雷达脉冲回波信号S′R(m,n)的幅度信息选取ISAR图像中目标所在的距离单元,即ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="661"/>得到目标回波信号S′R(l,n),其中,l表示目标所在的距离单元,l=i,i+1,...,j,i表示目标所在距离单元的起始位置,j表示目标所在距离单元的末尾位置。步骤4,利用多种准则从雷达脉冲回波信号S′R(m,n)中选择所需的距离单元。(4a)利用距离-多普勒算法对目标回波信号S′R(l,n)成像,得到ISAR图像I(l,n′),其中n′表示频率,n′=1,2,...,N,分别计算ISAR图像I(l,n′)每个距离单元的对比度c(l): c ( l ) = mean { [ | I ( l , n ′ ) | - mean [ | I ( l , n ′ ) | ] ] 2 } mean [ I ( l , n ′ ) ] , ]]>其中,|·|表示求模运算;(4b)将对比度c(l)进行降序排序,根据设定的对比度阈值μ1选择出c(l)的前L个值,从目标回波信号S′R(l,n)中提取这L个距离单元,作为初选距离单元,记为S′(l′,n),其中,l′是初选出来的L个距离单元,l′=1,2,...,L;(4c)计算初选距离单元S′(l′,n)的每个距离单元的对比度c(l′);(4d)分别计算初选距离单元S′(l′,n)每个距离单元的包络归一化方差v(l′),相位匹配二项式误差p(l′)以及每个距离单元的平均幅度a(l′): v ( l ′ ) = cov ( S ′ ( l ′ , n ) ) cov ( S ′ ( l ′ , n ) ) + [ mean ( S ′ ( l ′ , n ) ) ] 2 ]]> p ( l ′ ) = Σ n | angle ( S ′ ( l ′ , n ) ) - LSE ( angle ( S ′ ( l ′ , n ) ) ) | ]]> a ( l ′ ) = Σ n ( S ′ ( l ′ , n ) ) N , ]]>其中,cov(·)表示求方差运算,mean(·)表示求平均值运算,angle(·)表示求相位运算,LSE(·)表示最小二乘拟合直线;(4e)分别对包络归一化方差v(l′)进行升序排序,对相位匹配二项式误差p(l′)进行升序排序,对每个距离单元的平均幅度a(l′)进行降序排序,对每个距离单元的对比度c(l′)进行降序排序,依据排序后的v(l′),p(l′),a(l′)和c(l′)分别对初选距离单元S′(l′,n)中每个距离单元进行如下“投票”:(4e1)找到初选距离单元S′(l′,n)中每个距离单元在v(l′),p(l′),a(l′)和c(l′)中的排名顺序,分别对每个距离单元在v(l′),p(l′),a(l′)和c(l′)中排名顺序的名次进行相加,计算得到每个距离单元的名次总和;(4e2)对每个距离单元的名次总和重新进行升序排序,选择出前L′个值,从初选距离单元S′(l′,n)中选择出这L′个值对应的L′个距离单元,作为“投票”选择出的距离单元,记为S′(l′′,n),其中,l′′是新选择出来的L′个距离单元,l′′=1,2,...,L′,L′<L。步骤5,对“投票”选择出的距离单元S′(l′′,n)的每个距离单元进行局部多项式傅里叶变换lpft,得到调频率chirp(l′′)和每个调频率对应的归一化面积area(l′′):(5a)设定参数:α=-1,-0.999,-0.998,...,1,fn=16*N;其中,α表示调频率的搜索范围,fn表示傅里叶变换的点数;(5b)对“投票”选择出的距离单元S′(l′′,n)中第l′′个距离单元进行局部多项式傅里叶变换lpft(l′′,f,α),并计算第l′′个距离单元的图像对比度IC(l′′,α),并令l′′=1, IC ( l ′ ′ , α ) = mean { [ lpft ( l ′ ′ , f , α ) - mean [ lpft ( l ′ ′ , f , α ) ] ] 2 } mean [ lpft ( l ′ ′ , f , α ) ] , ]]>其中,lpft(l′′,f,α)为局部多项式傅里叶变换,计算公式为ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="700"/>该式可采用快速傅里叶变换实现,f表示频率,Tr表示脉冲重复周期,mean[lpft(l′′,f,α)]表示在频率维f上的求平均值运算;(5c)根据图像对比度IC(l′′,α)计算调频率chirp(l′′)和每个调频率对应的归一化面积area(l′′), chirp ( l ′ ′ ) = arg max α [ IC ( l ′ ′ , α ) ] ]]> area ( l ′ ′ ) = Σ α ( IC ( l ′ ′ , α ) ) , ]]>其中,ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="195"/>表示对表达式求最大时对应的α;(5d)更新l′′=l′′+1,转到步骤(5b),直至l′′=N。步骤6,对归一化面积area(l′′)进行升序排序,根据设定的归一化面积阈值μ2选择出area(l′′)的前K个值,从调频率chirp(l′′)中得到这K个值所对应的调频率chirp(k),从“投票”选择出的距离单元S′(l′′,n)中得到其对应的K个距离单元,从雷达回波信号S′R(m,n)中得到这K个距离单元对应的距离单元位置rangebin(k),其中,k表示最后选择出来的距离单元,k=1,2,...,K,K>1。步骤7,根据最小二乘方法,对调频率chirp(k)和其对应的距离单元rangebin(k)线性拟合,得到K个调频率估计的斜率ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="62"/>通过如下计算公式进行: a ^ k = Σ k = 1 K chirp ( k ) rangebin ( k ) - 1 K Σ k = 1 K chirp ( k ) Σ k = 1 K rangebin ( k ) Σ k = 1 K ( rangebin ( k ) ) 2 - 1 K ( Σ k = 1 K rangebin ( k ) ) 2 ]]>对调频率chirp(k-1)和其对应的距离单元rangebin(k-1)线性拟合,得到K-1个调频率估计的斜率ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="92"/>依次类推,直到对调频率chirp(2)和其对应的距离单元rangebin(2)线性拟合,得到2个调频率估计的斜率ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="59"/>最后计算所求斜率的平均值ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="46"/>即 a ^ = 1 K - 1 Σ k = 2 K a ^ k . ]]>步骤8,根据斜率平均值ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="44"/>得到垂直于成像平面的角速度ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="57"/> Ω ^ = c 2 f 0 | a ^ | ]]>其中,c为光速,f0为雷达射频。步骤9,依据角速度ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="57"/>对步骤(3a)中得到的ISAR图像I(l,n′)的方位维坐标n′乘上方位分辨率ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="136"/>即可完成ISAR图像的横向定标;再对ISAR图像I(l,n′)的距离维坐标l乘上距离分辨率ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="101"/>即可得到定标后的ISAR图像,表示为:ntent="drawing" img-format="tif" inline="yes" orientation="portrait" wi="409"/>其中,B为雷达信号带宽,Fr为脉冲重复频率,λ为波长。本发明的效果通过以下对实测数据的实验进一步说明:1.实验场景:本实验使用的实测数据是对雅克-42运输机观测的ISAR数据,机长实际长度为36.38m,雷达射频为5.52GHZ,脉冲重复频率为400HZ,读取数据间隔为12,读取的数据共有256个回波,每个回波包括256个采样点。2.实验内容:实验1,使用本发明方法对ISAR图像定标:2.1)读取雷达回波数据,并对该数据进行运动补偿;2.2)根据运动补偿后数据的幅度信息选取ISAR图像中目标所在的距离单元,得到目标回波信号;如图2所示为对运动补偿后的数据利用距离-多普勒算法成像的结果图。2.3)利用包络归一化方差,相位匹配二项式误差,每个距离单元的平均幅度,每个距离单元的频谱对比度等多种准则,从目标回波数据中选择所需的距离单元,如图3所示,图中圆圈标注的位置表示选择出来的距离单元。2.4)对选择出来的距离单元进行LPFT调频率估计,然后根据归一化面积准则选择最终的距离单元,如图3所示,图中星号标注的位置表示最终选择的距离单元。2.5)根据最小二乘方法,由最后选择出来的距离单元和其对应的调频率进行线性拟合,分别选择不同数目的距离单元和其对应的调频率,线性拟合得到多个直线斜率,最后求它们的平均值,作为最后得到的斜率。2.6)由得到的斜率估计垂直于成像平面的角速度,至此可对目标进行定标。实验2,本发明方法与常规方法的实验对比:本发明通过多种准则选择距离单元,然后进行LPFT调频率估计,在调频率估计之后对距离单元再次进行选择,最后通过最小二乘方法拟合多组数据,得到多个直线斜率,计算多个直线斜率的平均值作为最后得到的斜率,其选择的距离单元以及调频率的最小二乘拟合结果如图4所示。用常规方法仅使用包络归一化方差准则选择距离单元,直接进行LPFT调频率估计,得到距离单元和对应调频率的拟合斜率,其选择的距离单元以及调频率的最小二乘拟合结果如图5所示。3.实验结果分析:从图3中可以看出,初选出来的距离单元基本遍布目标所在的距离单元。由于以频域对比度为基础,靠近目标旋转中心,即目标重心能量较大的距离单元也被选择了出来,因为其靠近旋转中心,调频率较小,这部分距离单元不应该作为估计调频率的距离单元。星号标注的位置是最后选择出来的距离单元,其分布均匀,而且对靠近旋转中心的距离单元选择较少,能够较好的进行最小二乘估计。从图4和图5两幅图的对比中可以明显看出,在同样选择12个距离单元的情况下,常规方法通过最小二乘方法拟合直线的误差远大于本发明拟合直线的误差。通过两幅图拟合得到的斜率分别对目标定标,可以得到本发明对机长的定标结果为33.4619m,定标误差为8.02%,常规方法对机长的定标结果为30.1864m,定标误差为17.02%。通过使用本发明方法选择距离单元,定标结果能大幅改善。